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An analytic method for the solution of the evolution equation for GPDs is presented. The small 
x,£ asymptotics of GPDs are calculated. 
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One of the outstanding problems of quantum- 
chromodynamics is to understand how hadrons are built 
from quarks and gluons and how hadron properties re¬ 
sult from their microscopic quark-gluon structure. In 
recent years the formalism of generalized parton distri¬ 
butions (GPDs) was developed 00 , which offers the 
most economic and comprehensive possibility to describe 
the internal structure of hadrons. GPDs provide a gen¬ 
eral framework covering all relevant processes, from in¬ 
clusive ones (deep inelastic scattering) to exclusive ones, 
and, therefore, allow us to combine the information con¬ 
tained in different measurements in an optimal manner. 
In addition, GPDs link this information to fundamental 
elements of internal hadron structure, which cannot be 
directly determined by a measurement, like information 
on the total angular momentum of quarks or the spatial 
quark distribution in a fast moving hadron. 

For instance, knowledge of the transverse spatial struc¬ 
ture of quark and gluon distributions in a very high- 
energy proton is needed to describe proton-proton col¬ 
lisions at the LHC (CERN) which will start operating 
next year. The reason is that at LHC-energies several 
simultaneous hard interactions will occure in a single 
proton-proton collison. The rates for such events depend 
crucially on how quarks and gluons of different momen¬ 
tum fractions x are distributed in the transverse direc¬ 
tion. Without a quantitative understanding of this QCD- 
background e.g. the experimental search for physics be¬ 
yond the standard model will be less sensitive. 

In practice, the determination of GPDs is highly non¬ 
trivial, because they typically enter only in convolutions. 
Probably only global fits to all relevant experimental 0 
and lattice data 0 will be selective enough to really de¬ 
termine more than the most dominant GPDs. Such fits 
naturally rely heavily on Q 2 evolution. The evolution 
equations for GPDs are well known 0; however, analytic 
methods to solve them were lacking. As a result, one has 
to adhere to numerical methods which are usually cum¬ 
bersome and often unstable. Obviously, this whole field 
would profit enormously, if a faster and numerically more 
stable method to calculate the Q 2 evolution of GPDs 
could be found. In this contribution we present a simple 
algorithm, defined by Eqs. © - © which fulfills these 
demands for LO evolution. This algorithm is based on 


an analytic method for the determination of the GPD 
scale dependence, which relies heavily on the symmetries 
of the evolution equations. 

The analytic structure of GPD evolution is intrinsically 
related to that of Dokshitzer-Gribov-Lipatov-Altarelli- 
Parisi (DGLAP) evolution »G| for inclusive processes and 
Efremov-Radyushkin-Brodsky-Lepage (ERBL) evolution 
[j] for exclusive processes, because all three sets of evo¬ 
lution equations are only different facets of one and the 
same equation - the renormalization group equation for 
the corresponding composite operators. In this Letter 
we consider as concrete examples the gluon GPDs of the 
nucleon (H 9 (x,t;,t), E 9 (x,^,t)) related to the matrix el¬ 
ement of the twist-two gluon operator 

0(zi,z 2 ) = tr {G +9 (zin)G+(z 2 n)} , (1) 

but the method is general. We also neglect mixing 
with quark-antiquark operators. (For the exact defini¬ 
tions of the GPDs and kinematical variables, x, £, P + see 
Ref. 0.) Since all of these GPDs obey the same evolution 
equation we introduce the common notation, ip^(x). 

Because the derivation of our solution is rather techni¬ 
cal let us state already here the basic result. We find it 
crucial to treat the two kinematic regions |£| < |a;| (the 
DGLAP region) and \x\ < |£| (the ERBL region) differ¬ 
ently. Given some GPD parametrization at some input 
scale hi, one calculates the coefficients cj(j'), (j) from 

c?(j) = J ^ dxd lpA x ) d l K 1 (*!£!)]> ( 2 ) 

9,' _ i rVlil 

cf(j) = J^ dxdlpj{x)dl\ipf (±z|£|)] (3) 

and plugs these into 0 

00 / \ 
<P( 2 (x) = 0 (l£l - M)X c T ^ L ~ lb)p9 (^J 

Here L = and the anomalous dimension q(j) is 

specified later. The integration goes along the line par¬ 
allel to the imaginary axis such that 2 < Rej < 3. The 
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functions Pj(x), Qj(x) are expressed in terms of the Leg¬ 
endre functions of the first and second kinds Q 


pj(x) 


qj(x) 


(l-x 2 )Pj\(x), 


( 5 ) 

e‘-»Q7-Vx + )-e-»Q 7 _\(,-) 

2 i sin 7 jj 


where x± = x ± *0 D3- The sums and integrals in 0 
converge absolutely. As it is clear from 0 the first term 
contributes only to the ERBL region, while the integral 
gives contributions both to the DGLAP and ERBL re¬ 
gions. We notice also that function q 3 (x) = 0 for x < — 1, 
while for x > 1 and |x| < 1 it can be simplified to 


Qj(x) = (1 - ai 2 )Q J _ 2 1 (a;), x > 1 (7) 

<&(*) = (1 - x 2 )P~\(-x), \x\<l. (8) 

Z sin nj J 

It vanishes fast with j —> ±ioo. The formula 0 can 
be used both for numerical and analytic studies of GPD 
evolution. The main point is that if one chooses the GPD 
as a simple analytic form at the starting value /ii, the 
integrals in © and © can be done analytically. Then 
to restore GPD at the scale p 2 one has to evaluate the 
sum and the one-dimensional integral. We also stress 
that the form of Eq. 0 is completely determined by the 
symmetry properties of the evolution equation. 

Below we explain the main steps in the derivation of 
Eq. ©, the details will be given elsewhere. The evolu¬ 
tion equation for the GPD ip^(x) follows from that for the 
operator 0 . The corresponding kernel (see refs. Hi) 
is usually given in momentum space (i.e. as a function 
of x). For our purposes it is more convenient to use the 
coordinate-space formulation [lij 


/ OO 

dxe iP+xz ^(x) , (9) 

-OO 

<p(zi,z 2 ) = dte- iP+ ^+^<f>s(z 12 ), 


where Z\ and z 2 are real variables, and z\ 2 = z\ — z 2 . 
The LO evolution equation for z 2 ) reads 


d 

^ + 3i9 % 


(p{zi,z 2 ) = -N c — [H ip](zi,z 2 ). 
7r 


The integral operator H M can be cast into the form 


( 10 ) 


[H^] ( 2 : 1 , 22 ) = -4 da d(3(a/3 + 2aP)<p(zi 2 , z^) 


Jo Jo 
f 1 a 2 

+ / da — [2<p(zi, 2 2 ) - <p(*i 2 > z 2 ) 

Jo a 
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-^(2i,2fi)] + -^(2i, 2 2 ). (11) 

Here we use shorthand notations zf k = Zia + Zka, 
a = 1 — a. The Hamiltonian m commutes with the 


generators of collinear conformal transformations S a = 
(Si + S 2 ), a = ±,0; [H, S a ] = 0. The one particle 

operators, S%, are the generators of the SL( 2, R ) group 

S k = -<9fc, S£ = zldk + 2s k z k , S° = Zkd k + s k , (12) 

in the representation with spin Sk = 3/2. The corre¬ 
sponding finite transformations T(g ) maps the function 
f(z) of the real variable 2 onto another one according to 




1 

(cz + d) 3 


f 


( az + b\ 
\cz + d) ’ 


(13) 


where a, 6, c, d are the entries of the matrix g £ SL(2 , R). 
To solve the LO evolution equation this symmetry plays 
a crucial role. 

The GPD ip(zi, z 2 ) depends on two real variables 21 , z 2 
and transforms under SL( 2, R) transformations accord¬ 
ing to the tensor product of two representations d. 
Solving the evolution equation, one seeks for the com¬ 
plete system of functions that diagonalize the Hamilto¬ 
nian. Usually such a basis is provided by the eigenfunc¬ 
tions of the Casimir operator of the symmetry group. In 
the case under consideration the Casimir operator is 


JJ 2 = S+S~ + S°(S° - 1) = -z^did 2 zl 2 . (14) 


To solve the eigenvalue problem for the Casimir opera¬ 
tor o one should specify the scalar product. Our choice 
is restricted by the requirement of the SL( 2, R) invari¬ 
ance of the latter and by the condition for the “physical” 
GPD to be normalizable with respect to this scalar prod¬ 
uct. These requirements lead us to 


IMI 2 =/ dzidz 2 zf 2 |^( 2 i, 2 2 )| 2 = i J d^dx\d 2 ^(x)\ 2 . 

(15) 

One can check that d is invariant with respect to the 
transformations m and that the total spin operators 
S ±, ° are antihermitean. The Casimir operator d is 
self-adjoint and its eigenfunctions form a basis in the 
Hilbert space with the scalar product d- 

The eigenfunctions of the Casimir operator can be eas¬ 
ily found and they are in accordance with the decompo¬ 
sition of the tensor product of two representations d 
into irreducible ones Q The operator JJ 2 has both a 
discrete and continuous spectrum. The corresponding 
eigenfunctions can be expressed as 

PjL„ + i( 21 , 22 ) = VP^e~ iP+ ^+^^(z 12 ), (16) 
p ^ +ip (zi,z 2 ) = v / P+e- iP+ «(^+^)^’ ± (2 12 ), (1 7) 

where n £ Z+ and p £ R+. The functions 4')' ± (z) and 
4'/ ( 2 ) are defined for arbitrary complex j as follows 


4 ^( 2 ) 

4 - 5 ( 2 ) 


4-f(2±)-4-i_ J .(2±) 


(18) 


2 cos 7 rj L J 
e -zf 0-1/2) z -5/2j,_ i/2 (|P+£| 2 ) , (19) 
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where J u (z ) is the Bessel function and z± = ±z + iO. 
(Let us note that by construction \E'|’ ± (z) = (2).) 

The eigenvalue of the Casimir operator in both cases is 
E = j(j — 1). Taking the Fourier transform of (I Kill one 
finds that the eigenfunctions of the discrete and contin¬ 
uous spectrum have support in the ERBL and DGLAP 
regions, respectively. 

The decomposition of the arbitrary function ip(zi,z 2 ) 
into eigenfunctions of the Casimir operator thus reads 


ip{zi,z 2 ) 




{zi,z 2 ) 


( 20 ) 


pl/2+i 


' 1/2—i 


dj u c {j) af(j) Pf’ ± (z 1 ,z 2 ) 


where u>{j) = 2j — 1 and c o c (j) = (j — l/ 2 )cot 7 rj. The 
expansion coefficients a £ (j) and a £ (j) are given by the 
scalar product of ip(zi,z 2 ) with the eigenfunctions ll(il) 
and O- We get rid of the integral over £ and write 
down the expansion directly for the function $ 5 ( 2 ) © 
as follows 


$ 5 ( 2 ) 


^(2j-l)a £ (j)v|>|(2) 

3 =3 


( 21 ) 


1 

2 


/•5/2+ioo dj 
J 5 / 2 —ioo Sin 7 Tj 


{2j -l)af(j)^{z ± ). 


To derive 121 we made use of the explicit form of the 
functions m and the symmetry, a±(j) = a±(l—j), 
of the expansion coefficients. We notice also that the 
difference between the integrals in lt 20 l) and m can be 
calculated by residues and reproduces the missing terms 
(j = 1, 2) in the sum in the Eq. (1211) . 

Now we are ready to discuss the evolution of the func¬ 
tion $ 5 ( 2 ) driven by the Hamiltonian (Hill . The eigen¬ 
functions of the discrete spectrum of the operator J 2 are 
the eigenfunctions of the Hamiltonian with energies 


m = 2 


V’(j) - ^( 3 ) 


1 

j{j + 1 ) 


1 

U~ 1)(J - 2). 



Let us notice that E(j) is defined only for integer j > 3. 
Further, one can check that the eigenfunctions of the con¬ 
tinuous spectrum do not diagonalize H. This seems 
to contradict the commutativity of H and l 2 . The ex¬ 
planation is simple, the Hamiltonian 11 is not self-adjoint 
with respect to the scalar product CEJ> and thus does not 
have the same eigenfunctions as the operator JJ 2 E- 
Nevertheless, one can find functions which diagonalize 
the Hamiltonian (of course, they are not mutually or¬ 
thogonal and do not form a basis of the Hilbert space). 
Indeed, the equation l 2 ^(zi, z 2 ) = j(j— l)\H(zi, z 2 ) (here 
j is an arbitrary complex number) after separation of the 
factor e~ lP+ ^ Zl+z ' 2 \ turns into a second order differen¬ 
tial equation with the two independent solutions H/|( 2 ) 


and -(z). Because of the commutativity of the in¬ 
tegral operator HI and the differential operator l 2 , one 
concludes that 

Mj( 2 ) = A{j)**(z) + BU^tjiz). ( 22 ) 

Substituting <p{z\, z 2 ) = e~ iP+ £( Zl + S 2 )'kj( 2 ) into (TTT1) one 
finds that the integrals converges when Rej > 2. Next, 
let us notice that when the variables a and (3 run over the 
integration region the argument of the function varies 
from 0 to 2. Thus to fix the coefficients A(j) and B(j) it is 
sufficient to study the 2 —> 0 asymptotics of the r.h.s. and 
l.h.s. of d. In this case one can substitute T^z) by its 
leading term ~ z - 7-3 to get A(j) = E(j), B(j ) = 0 . Thus 
we have shown that the function e~ lP ^ Zl+ 22 )’F|( 2 ) di¬ 
agonalizes the Hamiltonian HI. 

Now we are in the position to solve the evolution 
equation for the GPD <I> £ (z). Indeed, the decomposi¬ 
tion m involves only the functions \I/|(z) which diag¬ 
onalize the Hamiltonian, and the integration follows the 
line Re j = 5/2. So we can apply the Hamiltonian H, 
Eq. Ill 111 , and change the order of integration between j 
and a, (3. Thus the solution of the evolution equation can 
be written in the form (cf. EUd) 


(*) = £>■-!) a T c?) L ~ 1{j) M (*) ( 23 ) 

3 =3 

2 J G sm nj « J 

Here 7 (j) = 2N c E(j)/bo , bo = ^N c — |n/. After Fourier 
transformation, Eq. d can be cast into the form 11- 

Let us discuss the properties of the solution ©. We 
represent GPD as <Pz(x) = <p\{ x ) + { p\ I { x )t where 
< p^(x ) is given by the sum in Eq. © and (x) by the 
integrals. The function <p^(x ) = 0 outside the ERBL 
region. The function ip^ix) = ipt(x) in the DGLAP 
region and does not vanish for |rc| < |£|. At the scale 
it can be calculated for |rr| < |£| in closed form, x ) = 
y> £ (£) + tp'^)(x 2 — £ 2 )/2£. This contribution is entirely 
due to the eigenfunctions of the discrete spectrum with 
j = 1,2 (see Eqs. (Pl> . d )- 

As it was stated earlier the sum and the integrals in 
© are absolutely convergent. Thus having calculated the 
coefficients c £ (j) and c £ ’ (j) at scale /i 1 one can easily 
restore the function at another scale /i 2 . 

Let us also note that the evolution of the coefficients 
c £ (j) and cjr (j) is not autonomous. To find the coeffi¬ 
cients at scale p 2 one should insert the function ip^ 2 (x ) 
into the Eqs. © and ©. Since the function qj{x) is 
not an eigenfunction of the Casimir operator one easily 
figures out, e.g. that the coefficient c| 2 (j) cannot be ex¬ 
pressed solely in terms of the coefficient c^{j). Rather, 








4 


after some algebra, one finds 


cfO'O = f 

J c 


dj 


1 


x L 


-7(7) 


<£ 1,+ (j)+ (-!)''^'"(j) . (24) 


c (f +3 - 1)0' - f) 




The basic Eq. 0 can be analyzed analytically in dif¬ 
ferent limits. As an example let us discuss the as ymp totic 
expansion of the GPD at large 112 (see Refs. a,!l3)- To 
this end we shift the integration contour from the line 
Re j = 5/2 to the line Re j = N — 1/2, N being an inte¬ 
ger. It follows from the Eqs. Q and © that for x > 1 the 
function qj(x) is analytic for Re j >2, while for \x\ < 1 
it has poles at integer j. Thus, shifting the contour one 
has to evaluate residues at these points. The integral 
over the new contour decreases faster than L~ 7( JV ~ 1 / 2 ). 
Keeping the first N — 3 terms from the sum and adding 
the contribution from the residues one finds 


N-l 

<&(*) = ©del - |*|) E (25) 

7=3 

Here Cj[_ 3 is the Gegenbauer polynomial and the coeffi¬ 
cient D^(j) reads 

/■via , /2 

D t 1 tt) = rj dxC j L 3 (x)ip^ 1 (x\^\), (26) 

where r 3 = 9/2(2 j - l)T(j - 2)/T(j + 2). 

As another illustration let us study the small £ behav¬ 
ior of GPDs. For small £ (we take £ > 0) the coefficient 
c+ ( j ) can be represented in the form 

cf’ + (i) = /C7)c^+(j)(2/£r 1 - (j - 1 - j) , (27) 

where /(j) = —-^r(j + l/ 2 )/P(j — 2 ) and the coeffi¬ 
cient c Ml,+ (j) = /p dxxi-^ip^L^x) + Lf(£ 2 ). Denoting 
by ip^ 2 ’ + (x) the integral corresponding to the term with 
a = + in the sum 0 and taking into account El) one 
finds for x £ 

<4 2 ’ + (x ) = j djx-i +2 c^'+{j)L-iW + <D(£ 2 ), 

(28) 

which corresponds to DGLAP evolution. It is interesting 
to compare the asymptotics in two regimes: £ —> 0 , x 
fixed and £ —> 0, x = a£. Using El and 0 one derives 
for the ratio of the asymptotics 


* + («> = % , 
^0 ( x ) 


M2,+ 


(X) e™f v {a + ) - e-™f v {a_) 


2 i sin7r^ 

1/2 


, (29) 


where a ± =a± * 0 , v = and 

}v{a) = e-^ ars(a) 2 Ei(^/2, v/2 + 1/2, ^ + 5/2|a" 2 ). 



a — */£ 


FIG. 1: The ratio A(a) = I? + (a) + R + {—a ) for L = 2.4 and 
2 ; = 10 -3 (upper curve), 10~ 4 , 10~ 5 and 10~ 6 (lowest curve). 


Taking into account the similar contribution from the 
a = — intergral in 0 one obtains for the total ratio 
R(a) = ip^ 2 a (x)/ipQ 2 (x) the result shown in Fig. ^ 

To summarize, we have developed an analytic method 
for solving the LO evolution equations for GPDs. This 
method is quite general and relies on the conformal sym¬ 
metry of the evolution equation. Our approach can 
also be used to elucidate the analytic structure of three- 
particle parton distributions which possess additional 
hidden symmetry 0 
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